# The following script was adopted from Thiele, Jan C., Kurth, Winfried and Grimm, Volker (2014) 
# 'Facilitating Parameter Estimation and Sensitivity Analysis of Agent-Based Models: 
# A Cookbook Using NetLogo and 'R'' Journal of Artificial Societies and Social Simulation 17 (3) 11 
# <http://jasss.soc.surrey.ac.uk/17/3/11.html>.


# a function to handle a simulation

simulate <- function(param.set, parameter.names, no.repeated.sim, nl.obj, trace.progress, iter.length, function.name) {

  if (length(param.set) != length(parameter.names))
  { stop("Wrong length of param.set!") }
  if (no.repeated.sim <= 0)
  { stop("Number of repetitions must be > 0!") }
  if (length(parameter.names) <= 0)
  { stop("Length of parameter.names must be > 0!") }

  # an empty list to save the simulation results
  eval.values <- NULL

  # repeated simulations 
  for (i in 1:no.repeated.sim)
  {
    # create a random-seed for NetLogo from R, based on min/max of NetLogo's random seed

    NLCommand("random-seed",runif(1,-2147483648,2147483647), nl.obj=nl.obj)
    
    NLCommand("setup", nl.obj=nl.obj)

    # set NetLogo parameters to current parameter values
    lapply(seq(1:length(parameter.names)), function(x) {NLCommand("set ",parameter.names[x], param.set[x], nl.obj=nl.obj)})

    # 12 months for start up
    
    NLDoCommand(1,"repeat 12 [go]", nl.obj=nl.obj)
    
    # 144 months for 12 years model simulation (12 x 12)
    
    cal.crit <- NLDoReport(12,"repeat 12 [go]",c("unemployment-rate","ltu-rate", "theta"), 
                           as.data.frame=T, df.col.names=c("u","ltu", "theta"), nl.obj=nl.obj)


    # calculate calibration criteria

    umean.criterion <- mean(cal.crit$u) 
    theta.criterion <- ( sd(cal.crit$theta) / mean(cal.crit$theta) )
    ltu.criterion <- mean(cal.crit$ltu)


    # merge calibration criteria
    calibration.criteria <- c(umean.criterion, theta.criterion, ltu.criterion)

    # append to former results
    eval.values <- rbind(eval.values,calibration.criteria)
  }

  # print the progress if requested
  if (trace.progress == TRUE)
  {
    already.processed <- get("already.processed",env=globalenv()) + 1
    assign("already.processed", already.processed, env=globalenv())
    print(paste("processed (",function.name,"): ", already.processed / iter.length * 100, "%", sep = ""))
  }

  # return the mean of the repeated simulation results
  if (no.repeated.sim > 1)
  {
    return(colMeans(eval.values))
  }
  else {
    return(eval.values)
  }
}
